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I.  Introduction 

The  United  States  Navy  uses  S5nthetic  coatings  for  various  purposes 
on  a  variety  of  ships.  These  coatings  are  typically  viscoelastic  polymers  that 
are  prepared  and  applied  to  the  vessels  in  the  shipyard  during  some  phase 
of  construction. 

The  coating  begins  with  its  constituent  parts  separated,  some  of 
which  are  solid  and  some  liquid.  The  components  are  mixed  and  the 
product  applied  to  the  hull  is  a  highly  viscous  liquid.  A  crawling  apparatus 
moves  from  top  to  bottom  depositing  a  strip  of  coating  to  the  hull.  The 
crawler  then  moves  over  and  applies  a  strip  adjacent  to  the  previous  one. 
The  liquid  cures  very  quickly  into  a  solid  form  thus  retaining  the  desired 
thickness  on  the  hull.  The  proportions,  mixing,  and  application  of  the 
constituents  must  be  done  with  a  high  degree  of  accuracy  and  consistency 
in  order  for  the  final  coating  layer  to  meet  performance  expectations.  Any 
errors  in  the  mixing  or  delivery  process  will  result  in  affected  material 
properties  such  as  Young's  and  shear  modulus. 

Currently,  a  coring  method  is  used  to  determine  the  quality  of  the 

pour.  In  this  method  a  sample  of  the  coating  is  removed  from  the  layer 

after  the  pour  is  complete  and  is  then  sent  to  a  laboratory  where  various 

material  property  measurements  are  made  on  the  sample.  It  is  assumed 

that  the  properties  of  this  core  are  representative  of  the  layer  as  a  whole. 

The  area  from  where  the  sample  was  removed  is  then  patched  or  filled 

with  leftover  material.  There  are  several  disadvantages  to  this 
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methodology:  1)  it  is  destructive  2)  it  is  a  local  or  point  measurement  only 
3)  there  is  a  delay  caused  by  the  need  to  send  the  sample  to  a  laboratory. 
The  first  disadvantage,  the  destructive  nature  of  the  methodology,  is 
undesirable  because  it  creates  a  discontinuity  in  the  surface  and  introduces 
additional  interfaces  in  the  coating  layer.  The  second  disadvantage,  the 
local  nature  of  the  measurement,  is  undesirable  because  it  returns  the 
properties  for  a  very  small  region  of  the  coating.  If  the  properties  varied 
as  the  pour  progressed,  the  properties  far  away  from  the  core  location 
might  be  significantly  different.  The  third  disadvantage,  delay,  is 
undesirable  because  it  does  not  provide  for  quality  assurance  during  the 
process.  That  is,  if  the  process  suddenly  began  to  stray  outside  of  desired 
parameters,  the  operators  may  not  know  until  days  after  the  pour  is 
completed.  Only  later  would  the  failure  be  known  and  a  large  balance  of 
the  pour  would  fail  acceptance  criteria. 

This  research  proposes  a  different  method  for  determining  material 

properties  that  avoids  the  aforementioned  disadvantages  of  the  current 

system.  The  proposed  system  would  avoid  the  destructive  disadvantage 

by  using  a  non-contact  sensor  for  detection  of  signals  generated  by  a 

surface  contacting,  but  not  penetrating,  shaker.  The  proposed  system 

would  avoid  the  local  disadvantage  by  making  measurements  that  would 

test  material  properties  over  a  significant  spatial  range  and  repeating  these 

measurements  at  many  locations  over  the  entire  surface.  The  proposed 

method  would  avoid  the  delay  disadvantage  by  processing  the  data 
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quickly  enough  to  alert  the  operators  in  time  to  halt  and  fix  the  process 
very  shortly  after  the  pouring  system  strayed. 

The  work  presented  in  this  report  represents  progress  made  up  to 
the  point  of  termination  of  the  project  due  to  funding  considerations. 
Although  only  43%  of  the  contract  value  was  funded,  the  major  goal  of 
identifying  a  method  for  determining  material  properties  has  accomplished 
and  application  of  the  technique  developed  here  to  several  other  related 
problems  has  been  identified  and  one  of  those  investigated,  (see  "VII. 
Comments  and  Conclusions"  for  a  listing  of  these  additional  applications 
and  "VI.  Results"  for  a  presentation  of  the  data  taken  and  modelling  done 
in  pursuit  of  one  of  these  alternate  applications).  These  alternate 
applications  were  discovered  and  the  one  pursued  without  a  prepared 
out-of-specification  sample  which  was  one  of  the  original  contract  items 
not  delivered  due  to  funding  considerations. 
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II.  Optimization  Approach 

The  technique  proposed  here  for  the  determination  of  material 
properties  is  the  use  of  an  optimization  process.  In  this  process  measured 
data  is  quantitatively  compared  to  the  results  from  a  theoretical  model 
whose  inputs  include  the  properties  of  the  subject  material.  An  iterative 
process  is  used  to  intelligently  vary  the  input  parameters  in  order  to 
improve  the  correlation  of  the  theoretical  model  results  to  the  measured 
data,  thereby  determining  the  material  properties. 

The  object  of  comparison  is  a  transient  sinusoidal  signal.  Although 
the  theoretical  model  is  better  suited  to  a  continuous  wave,  a  transient 
was  necessary  due  to  the  finite  dimensions  of  the  test  sample.  The 
technique  calls  for  a  single  cycle  of  the  pulse  to  be  delivered  to  the  free 
surface  of  the  sample  and  its  propagation  recorded  by  measuring  the 
surface  velocity  at  discrete  (small)  points  along  the  surface  of  the  sample 
radiating  outwards  from  the  epicenter.  The  measurement  system  records 
data  in  the  time  domain  and  then  the  data  processing  transforms  the 
recorded  data  into  the  frequency  domain. 

The  theoretical  model  takes  the  geometry  of  the  sample,  its 
constraints,  the  frequency  of  oscillation,  position  of  interest,  and  the 
material  composition  as  its  input  variables.  The  output  is  the  surface 
velocity  at  the  position  of  interest  and  the  inputted  frequency.  This  result 
is  a  continuous  wave  result,  so  a  convolution  with  the  input  signal  is 

necessary  to  predict  the  transient  response. 
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At  this  point  in  the  procedure  there  are  experimental  and  theoretical 
data  for  each  position  and  frequency  in  the  range  of  interest  which  is 
determined  by  the  power  spectrum  of  the  input  signal.  The  experimental 
data  are  pairs  of  complex  velocities,  one  being  the  component  normal  to 
the  surface  of  the  polymer  layer  and  the  other  being  parallel  to  the 
surface.  An  automated  optimization  routine  then  varies  the  material 
properties  and  searches  for  a  smaller  error  between  the  experimental  and 
theoretical  results  until  a  predetermined  convergence  value  is  achieved.  A 
schematic  of  the  optimization  routine  is  shown  in  figure  1. 


[time  of 
flight]  ■ 


approx.  cL 


[Ray  Model 
+Simplex] 


better  cL* 


[End] 

figure  1;  optimization  routine  schematic 
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III.  Experimental  Setup 

For  this  research  the  sponsor  provided  a  test  specimen  of  the  subject 
polymer.  It  was  a  two  inch  thick  slab  with  thirty  inch  sides.  The  slab  was 
attached  using  epoxy  to  a  two  inch  thick  steel  substrate  plate  with  thirty- 
six  inch  sides. 

(see  figure  2) 


figure  2:  sample  moimted  on  steel  backing  plate  with  shaker  and  LDV 
probe  positioned  on  the  free  surface  of  the  test  sample 

The  detection  system  was  a  Laser  Doppler  Velocimeter  with  a  two- 
axis  detection  head  capable  of  measuring  the  in-plane  (parallel)  velocity 
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and  out-of-plane  (normal)  velocity  in  near  simultaneous  fashion,  (see 
Appendix  A  for  details  of  the  LDV  system).  The  head  used  fiber  optic 
connections  and  was  moimted  on  a  computer-controlled  positioning  device 
that  allowed  for  precise  spacing  between  the  measurement  points  and 
convenient  focusing.  Time  domain  waveforms  were  digitized  using  a 
Tektronix  2430A  oscilloscope  and  a  GPIB  connection  used  to  transfer  the 
data  to  a  computer  for  storage  and  analysis. 

A  Bruel  &  Kjaer  model  4810  mini  excitation  shaker  was  used  to 
generate  the  waves  input  to  the  sample.  A  Wavetek  275  Signal  generator 
created  the  electrical  signal  that  drove  the  shaker  using  one  cycle  of  a  sine 
wave  at  2  Kilohertz.  An  accelerometer  was  mounted  to  the  head  of  the 
shaker  and  served  as  the  contact  surface  for  imparting  the  disturbance 
into  the  slab.  See  figure  3.  The  accelerometer  signal  was  recorded  and 
digitized  for  use  during  processing. 
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B&Kmini 
excitation  shaker 


figure  3:  illustration  of  the  B&K  shaker  resting  on  the 
coating  slab  surface 

The  shaker  was  placed  approximately  in  the  middle  of  the  slab  and 
measurement  locations  identified  radially  extending  from  the  epicenter,  as 
shown  in  figure  4.  Measurement  point  1  is  0.064  meters  from  the  epicenter 
and  subsequent  locations  are  spaced  at  0.5"  (0.0127  m)  intervals.  Accurate 
spacing  was  achieved  by  using  the  stepper  motor  controls  attached  to  the 
positioning  arm  to  which  the  probe  head  was  attached. 
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IV.  Theoretical  Modelling 

Several  theoretical  models  were  investigated  and  used  in  varying 
degrees  for  this  research.  A  simple  ray  acoustics  model  was  used  for  basic 
investigation  and  seed  values  for  the  more  complex  model  (see  V. 
Optimization  Process).  A  model  by  Miller  and  Pursey^  was  used  for  initial 
comparisons  and  understanding  of  the  basic  wave  structure.  The  limitation 
of  this  model  is  that  it  considers  only  the  case  of  a  free  half-space.  For  the 
final  processing  a  model  developed  by  Weaver^  was  used  since  it  was 
completely  applicable  in  terms  of  geometry  and  input  variables. 

The  geometric  ray  acoustics  model  used  a  theoretical  input  signal 
that  closely  approximated  the  actual  input.  One  cycle  of  a  two  kilohertz 
sine  wave  was  generated  in  the  Matlab®  computational  environment.  To 
improve  its  relation  to  the  actual  signal  a  Hanning  window  was  applied, 
effectively  removing  the  first-derivative  discontinuity  at  the  beginning  and 
end  of  the  cycle.  This  discontinuity  does  not  exist  in  the  physical  signal 
since  the  acceleration  of  the  shaker  head  is  not  instantaneous.  This 
waveform  was  then  replicated  three  times,  each  using  a  different  delay 
and  decay  based  on  the  appropriate  wave  speed  and  attenuation  value  (in 
addition  to  spherical  spreading).  The  three  formed  cycles  correspond  to 
the  direct  arrival  of  the  compressional  wave,  the  direct  arrival  of  the  shear 
wave,  and  the  arrival  of  the  compressional  wave  after  reflection  off  of  the 

^  G.F.  Miller  and  H.  Pursey,  "The  field  and  radiation  impedance  of  mechanical  radiators  on  the  free 
surface  of  a  semi-infinite  isotropic  solid".  The  National  Physical  Laboratory,  Teddington,  Middlesex. 

^  "Transient  ultrasonic  waves  in  a  viscoelastic  plate:  Theory".  Richard  Weaver,  Wolfgang  Sachse,  and 
Lin  Niu,  Department  of  Theoretical  and  Applied  Mechanics,  Cornell  University,  Ithaca,  New  York 
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coating-steel  interface.  The  polarity  of  the  formed  signals  was  matched  to 
the  directionality  of  the  LDV  probe,  (see  Appendix  A  for  details  of  the 
LDV  system).  The  formed  cycles  were  then  combined  to  create  a 
composite  signal. 

The  Miller  and  Pursey  model  is  an  adaptation  of  Lamb's  method 
designed  to  furnish  definite  integral  representations  of  the  field  at  an 
arbitrary  point  in  a  semi-infinite  isotropic  solid  due  to  prescribed  periodic 
stresses  on  the  free  surface.  The  equation  of  wave  propagation  in  an 
isotropic  solid  is  given  here: 

—  —  u 

CijVV.u  -  C44VX Vx  u  = 

dt 

equation  1:  wave  propagation  in  isotropic  solid 
where 

u  -  displacement 

cll  -  compressional  elastic  constant 
c44  -  shear  elastic  constant 
p  -  density  of  the  medium 

Solution  of  this  equation  solved  for  the  boimdary  condition  of  a 
circular  disk  vibrating  normally  to  the  free  surface  yields  the  following 
equations. 
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'  <^44  i  ^o(?)  I 


equation  2:  Miller  and  Pursey  result  for  circular  disk 
vibrating  normally  to  the  free  surface  (out  of  plane  displacement) 

where 


a  -  radius  of  driving  piston 
^  -  integration  variable 


^2 


2(1 -<t) 

1-2ct 


where 

G  =  Poisson' s  ratio  for  the  medium 


and  F„  =  (2f  ^  -  /i")'  -  4?^  ^  -  l)^(f  ^  ) 

The  equation  for  the  radial  displacement  is  similar.  Note  that  the 
input  variables  are  the  complex  moduli,  driving  piston  radius,  and  the 
radial  and  axial  position  of  the  evaluation  point.  An  example  of  the  results 
from  this  model  is  shown  in  figure  5. 
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Miller+Pursey  model,  r=0.0894m 


figure  5:  result  of  the  Miller  and  Pursey  model 
for  r=0.0894  m,  in  plane  and  out  of  plane  displacement 


This  model  was  used  to  gain  an  understanding  of  how  the  changes 
in  modulus  value  affected  the  shape  of  the  output  signal.  This  was  useful 
in  examining  the  structure  of  the  measured  waveforms  in  order  to  develop 
analysis  procedures,  particularly  the  use  of  the  ray  model  in  determining 
approximate  seed  values  for  the  Simplex  optimization  procedure,  (see  V. 
Optimization  Process). 
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The  third  model  used  was  the  one  developed  by  Weaver,  Sachse, 
and  Niu  with  the  stated  goal  of  determining  the  nearfield  response  of  a 
thick  viscoelastic  plate  at  epicenter,  off-epicenter,  and  on  the  same  side  of 
the  plate  to  a  point  step  load  acting  normal  to  the  surface.^  In  their 
derivation,  the  test  "object"  is  a  homogeneous  infinite  flat  slab  of  thickness 
2h  (=H)  with  free  surfaces  at  z=±h  subjected  to  a  concentrated  downward 
force: 

r^=~Femr) 


where 

0{t)  =  step  function  in  time 
5(r)  =  delta  function  in  space 

equation  3:  force  applied  to  free  surface  in  the  Weaver  model 


The  displacement  potentials  (j)  and  If/  are  defined  by 

dr  dz 


30 

w  =  ^  + 


1  d(ry/) 


dz  r  dr 

equation  4:  displacement  potentials  for  Weaver  model 

where 

u  =  radial  displacement  field 


^  "Transient  ultrasonic  waves  in  a  viscoelastic  plate:  Theory".  Richard  Weaver,  Wolfgang  Sachse,  and 
Lin  Niu,  Department  of  Theoretical  and  Applied  Mechanics,  Cornell  University,  Ithaca,  New  York 
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w  =  vertical  displacement  field 


The  solution  for  the  vertical  displacement  field  is  given  in  the 
frequency  domain  in  terms  of  the  Fourier  transformed  normal 
displacement. 


’'s  0 


_sin[ph)sm{ah)  ^  cos[ph)cos{ah) 


R. 


R. 


dk 


equation  5:  result  of  inverse  Hankel  transform  for  the  vertical 

displacement  field 


This  equation  is  then  integrated  over  all  values  of  k  to  achieve  the 
final  answer.  There  is  some  difficulty,  however,  in  that  the  high 
wavenumber  portion  of  the  integral  does  not  necessarily  converge 
absolutely.  To  overcome  this  convergence  problem  an  analytical  expression 
is  used  to  cover  this  region  of  the  integration. 


equation  6:  analytical  expression  for  high  wavenumber 
portion  of  the  integrand  in  Weaver's  model 

To  evaluate  this  solution,  the  appropriate  input  was  needed  in 

addition  to  material  properties  and  geometries.  Additionally,  a  boundary 

condition  had  to  be  changed  from  the  free-free  condition  previously 

existing  to  a  rigid  boimdary  condition  assumed  on  the  lower  surface 
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where  the  slab  was  attached  to  the  steel  plate. 

The  theoretical  input  used  a  step  function  in  time  which  was 
inappropriate  for  the  experimental  conditions  which  used  a  one  cycle  sine 
wave.  To  utilize  a  more  accurate  input  force  for  the  equations,  the 
recorded  and  digitized  input  signal  was  used  in  place  of  the  step  fimction. 
Additionally  the  boimdary  condition  at  the  lower  surface  of  the  slab  was 
not  appropriate  for  the  experiment  that  was  conducted.  The  theory  treats 
the  lower  surface  as  free  whereas  the  experimental  conditions  have  it  as 
(nearly)  rigid.  When  this  change  of  boundary  conditions  is  addressed  and 
the  theory  reformulated,  the  integrals  to  be  determined  are  as  follows. 


u,{r,co) 

where 


l(2^^-kl)d,-2^^pd 


dj  —  —CC  + 


2a 


,2  , 2'  cos(a'L)  cos(jS'L)  + 


2a 


2f  -^2 


2f -^2 


ap  sin(cif 'L)  sin(  jS'L) 


2a 


2^' 


-  [aj8cos(/?'L)  cos( or'L)  +  sin(/?'L)sin(«'L) 


a  =  -k^  =  ia' 

^  is  the  Hankel  transform  variable 
a  is  the  drive  piston  radius 

Similarly,  the  vertical  displacement  integral  is  given  next. 
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.2  la 


M,(r,<0)  =  J 


a 


•^iK) 


0(2^^ -k^)- 213^^ 


2a 


l^-kl 


equations  7:  final  Weaver  integrals 


The  integrals  were  then  evaluated  over  a  range  of  frequencies  for 
which  the  measured  input  signal  had  energy  (spectral  range).  These  result 
were  then  inversely  transformed  into  the  time  domain  for  qualitative 
comparison. 

A  sample  time  domain  result  is  shown  in  figure  6. 
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am|Gii1ude 


Weaver  mod  el  result,  vertical  daplacement 


figure  6:  Weaver  model  result  for  vertical  displacement 
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V.  Optimization  Process 

In  order  to  systematically  vary  the  input  parameters  to  improve 
correlation  between  measured  and  calculated  results  a  Simplex  routine  was 
implemented.  This  algorithm  operates  by  changing  the  size  and  shape  of 
an  object  in  an  n-dimensional  solution  space  where  n-1  equals  the  number 
of  variables  being  solved  for.  The  process  requires  an  initial  n-number  of 
solution  sets  as  a  starting  point  but  the  resulting  shape  need  not  boimd 
the  actual  solution  thereby  making  it  a  very  robust  optimization  process. 
The  routine  then  computes  a  single  performance  metric  for  the  current 
solution  sets.  If  the  metric  does  not  meet  a  prescribed  success  tolerance,  it 
automatically  reshapes  the  solution  volume  by  moving  the  vertex  which 
had  the  poorest  performance  in  the  previous  iteration.  Eventually,  as  the 
solution  volume  converges  around  the  solution,  the  solution  volume 
shrinks  and  the  vertices  converge  to  the  final  value.  For  more  information 
on  this  technique,  refer  to  Numerical  Recipes:  The  Art  of  Scientific 
Computing.  (New  York:  Cambridge  University  Press,  1986).  Appendix  C 
of  this  document  lists  the  code  used  and  a  sample  Simplex  solution. 

This  technique  was  chosen  for  its  robust  nature  and  its  ease  of 

programming.  Test  cases  were  performed  on  two-dimensional  problems  so 

that  the  progress  of  the  algorithm  could  be  monitored  by  plotting  the 

solution  volume,  in  this  case  a  triangle,  on  an  x-y  plot.  It  was  then  easy  to 

build  confidence  in  the  routine  as  its  convergence  was  monitored  for  this 

test  case.  A  disadvantage  of  the  routine  is  that  it  requires  an  initial  solution 

24 


volume  near  but  not  necessarily  encompassing  the  actual  solution  in  order 
to  avoid  becoming  trapped  in  a  local  minima. 

To  initiate  the  Simplex  routine  for  the  case  at  hand  a  procedure  was 
needed  to  determine  approximate  material  parameters  that  does  not 
require  an  initial  guess.  To  do  this,  a  simple  geometric  acoustics  ray  model 
was  employed  in  conjunction  with  time  of  flight  measurements. 

The  Ray  model  is  a  simplification  that  does  not  include  all  of  the 
physics  involved  in  the  situation.  However,  it  is  sufficient  as  a  means  of 
calculating  initial  guess  values  for  the  Simplex  procedure  operating  on  the 
Weaver  model.  The  Weaver  model  (in  conjimction  with  a  variation  on  it 
for  a  half-space  by  Miller  and  Pursey)  showed  that  the  fastest  wave 
traveling  along  the  free  surface  of  the  slab  is  directly  related  to  the  plane 
wave  modulus.  Given  this,  it  is  reasonable  to  use  the  first  half-cycle  of  the 
recorded  signal  to  determine  the  longitudinal  speed  of  the  material.  By 
measuring  the  result  at  five  locations  extending  along  a  radial  away  from 
the  epicenter  it  is  possible  to  also  get  an  estimate  of  the  attenuation.  To  do 
this,  the  ratio  of  the  amplitudes  of  the  successively  recorded  waves  were 
calculated  and  the  decrease  due  to  spherical  spreading  removed.  The 
remaining  decrease  was  attributed  to  viscoelastic  loss  in  the  material.  By 
doing  this  calculation  in  the  frequency  domain  on  a  per-frequency  basis  the 
effects  of  dispersion  are  automatically  accoimted  for  since  the  results,  both 
speed  and  attenuation,  are  at  discrete  frequency  points. 

The  wave  constituents  considered  in  the  geometric  ray  model  were 
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the  direct  longitudinal  wave,  the  direct  shear  wave,  and  the  reflection  of 
the  longitudinal  wave  off  of  the  interface  with  the  steel  plate  at  the  bottom 
of  the  slab.  A  sample  result  is  shown  in  figure  7. 
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figure  7:  Ray  model  result 


For  comparison,  the  ray  model  was  run  and  its  results  compared  to 
the  LDV  data  taken  for  the  third  measurement  point,  shown  in  figure  8. 
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figure  8:  ray  model  result  overlaid  with  LDV  for  point  3 
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VI.  Results 

The  shaker  used  to  excite  the  coating  slab  had  an  accelerometer 
mounted  to  its  head  and  was  therefore  in  direct  contact  with  the  coating 
during  the  excitation.  The  recorded  signal  from  the  accelerometer  was 
used  in  the  theoretical  calculations  and  is  shown  in  figure  9. 

The  LDV  system  was  used  to  record  waveforms  at  the  measurement 
points  for  both  in-plane  and  out-of-plane  motion.  Those  recorded  signals 
are  shown  in  figures  10  and  11. 


Accelerometer  Drive  Signal 
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In  Plane  motion  -  LDV 


0  0.5  1  1.5  2  2.5  3  3.5  4 

time,  seconds  xIO’^ 

figure  10:  In  Plane  motion  signals  detected  by  LDV 

In  figure  10  there  are  five  plots  shown  for  data  taken  in  the  time 
domain.  The  first  location  is  at  a  position  located  0.064  meters  from  the 
epicenter,  (refer  to  figure  4,  page  13).  Each  subsequent  point  (2-5)  is  taken 
at  a  location  0.0127  m  further  away  from  the  epicenter  on  a  line  extending 
radially  away  from  the  epicenter. 
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Out  of  Plane  motion  -  LDV 


figure  11:  out-of-plane  motion  signals  detected  by  LDV 


Weaver  model  in-plane  and  out-of-plane  results  for  all  five  location 
points  are  presented  in  figure  12. 

The  results  of  the  Simplex  convergence  are  shown  in  figure  13  and 
14. 
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R  Hank  el  result 


Z  Hankd  result 


figure  12:  Weaver  results  for  R  (in-plane)  and  Z  (out  of  plane) 
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figure  13:  Simplex  convergence  for  Longitudinal  speed 
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figure  14:  Simplex  convergence  for  Shear  speed 


For  comparison,  shear  modulus  data  was  obtained  from  an 
experiment  using  a  different  method  of  measuring  material  properties  on 
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the  same  subject  material/  A  comparison  table  is  presented  below. 


Hankel  Model  +  Simplex 

Conventional  Method 

frequency,  Hz 

977 

1100 

Shear  speed,  m/s 

153  +i  21.7 

126  +i  29.4 

table  showing  a  comparison  of  the  results  for  the  shear  speed 
from  a  conventional  method  to  the  one  presented  here 


Please  see  "Comments  and  Conclusions",  section  VII,  for  a  discussion 
of  these  results  and  their  comparison. 

Although  the  original  goal  of  the  program  was  to  develop  a 
technique  for  measuring  material  properties,  a  tangential  application  of  the 
experimental  technique  developed  here  was  briefly  explored.There  is  great 
interest  in  the  detection  of  voids  or  occlusions  in  a  poured  layer  which  can 
come  from  either  air  bubbles  or  coagulated  deposits  of  a  constituent 
material.  In  either  case  the  acoustical  transmission  through  the  coating 
layer  is  adversely  affected  by  the  presence  of  occlusions.  Due  to  the 
difference  in  density  and  elastic  properties,  sound  waves  will  scatter  off  of 
occlusions  and  thus  provides  a  means  of  detection.  The  LDV  system  was 
used  to  measure  wave  propagation  in  a  pair  of  samples  made  from  the 
same  material  used  in  the  body  of  this  research.  A  0.5"  wide  strip  of  the 

^  Willis,  Richard  Lance.  "Non-Invasive  Characterization  of  Micro>voided  Polymers  Under  Controlled 
Static  Pressure  and  Temperature  Using  Laser  Doppler  Vibrometry".  Presented  as  a  PhD  thesis  at  the 
Georgia  Institute  of  Technology,  School  of  Mechanical  Engineering,  December  1999. 
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material  was  cut  off  of  the  main  slab  and  mounted  on  a  steel  block  with 
the  same  thickness  as  the  moimt  plate  for  the  slab.  A  hole  was  then  drilled 
through  the  strip  to  simulate  an  occlusion  being  present  in  the  material. 
The  wave  propagation  through  the  region  containing  the  occlusion  was 
compared  to  propagation  through  unaltered  material.  The  results  show  a 
noticeable  difference  and  therefore  an  ability  to  detect  their  presence.  See 
figures  15-17. 

Although  this  application  was  not  explored  further  it  does  show  the 
versatility  and  potential  usefulness  of  the  technique.  In  addition  to  the 
occlusion  detection,  ideas  for  additional  applications  include  bonding 
integrity  and  seam  integrity  analysis. 


measurement 
locations 


steel  base 


point  force  input 
(transient) 


"Sdid 

side 


"Cavity" 
side 


material  strip 
(15‘^2‘'x0.75") 


0.25"  hole 
drilled  through 
at  a  depti  equal  to 
1/3  of  strip  height 


figure  15:  occlusion  detection  test  sample 
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12  3  4 


"Solid  side" 


'Deep"  cavtty 


"Shatlou/’  cavity 


figure  16:  nomenclature  for  occlusion  detection 
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figure  17:  results  for  shallow  cavity  detection 

Although  no  analysis  will  be  presented  here,  the  presence  of  the 
occlusion  can  be  detected  qualitatively  by  observing  the  changes  in  the 
front  end  of  the  recorded  waveforms  from  the  "Solid  Side"  to  the  "Cavity 
Side".  This  difference  is  also  plotted  in  the  third  plot,  "Subtraction". 
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VII.  Comments  and  Conclusions 

The  data  presented  in  the  results  section  reflects  the  best  results 
obtained  at  the  time  of  the  suspension  of  research.  Due  to  the  nature  of 
the  experimental  work  it  was  not  possible  to  test  a  known  material  and 
thereby  confirm  the  accuracy  of  the  technique  presented.  Likewise,  the 
usefulness  of  the  comparison  to  the  conventional  method  is  somewhat 
limited  in  that  it  was  not  tested  against  a  known  material  either. 

This  technique  was  designed  to  be  used  in  a  demanding  situation 
(bonded  on  one  surface)  and  the  success  achieved  is  encouraging.  There 
are  three  main  items  which  should  be  addressed  if  better  accuracy  is 
desired  from  this  technique  in  the  future. 

First,  measurements  should  be  made  using  a  force  gauge  in  the 
delivery  chain.  The  Weaver  model  uses  a  stress  input  and  the  results 
obtained  in  this  work  used  the  acceleration  signal  for  input.  Although 
closely  related,  the  differences  in  these  measured  quantities  do  create 
error. 

Second,  a  better  understanding  of  the  Simplex  method  should  be 
acquired.  For  complicated  solution  spaces  like  the  one  present  in  the 
integral  solutions  here,  the  initial  input  formulation  is  important  so  as  to 
avoid  false  minima.  Some  investigation  should  be  made  into  the  merits  of 
either  a  successive  two-parameter  search  or  a  singular  four-paramter 
search.  Again,  these  two  methods  may  have  different  behaviors  in 

searching  the  complicated  solution  space. 
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Third,  some  investigation  into  the  subtleties  of  the  numerical 
integration  should  be  made.  The  integrands  in  the  present  work  exist  in 
the  complex  plane  with  complicated  behaviors.  The  numerical  integration 
interval  must  be  properly  chosen  and  it  is  possible  that  a  non-regular 
stepping  would  be  advisable  in  order  to  properly  resolve  areas  in  which 
integrand  values  change  rapidly.  Additionally,  the  transition  from  the 
standard  integral  to  the  asymptotic  form  must  be  done  at  an  optimum 
point.  Further  investigation  of  these  issues  may  result  in  better  output 
from  the  theoretical  model. 

Despite  these  areas  of  possible  weakness,  the  combination  of 
experimental  and  theoretical  work  presented  here  is  considered  successful. 
The  Laser  Doppler  Vibrometry  system  worked  well  in  obtaining  the 
desired  data  for  the  technique.  The  theoretical  models  produced  results 
that  are  expected  and  logically  consistent.  The  output  from  the  models 
compares  well  to  the  experimental  data  on  a  qualitative  level,  particularly 
shape  and  features  of  the  waveforms.  The  quantitative  results  compared 
reasonably  well  with  the  other  method  and  fit  into  the  library  of  results 
for  this  type  of  material  quite  well. 

Experience  in  fine-tuning  design  and  operation  of  Laser  Doppler 

Vibrometer  systems  was  gained.  As  the  system  used  for  this  work  has 

also  been  used  for  previous  and  subsequent  experimental  work,  the 

hardware  constantly  undergoes  improvements.  This  research  effort  raised 

the  level  of  sophistication  and  furthered  knowledge  of  LDV  hardware 
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design  and  usage. 

Development  of  the  theoretical  models  required  learning  a 
considerable  amount  about  wave  propagation  in  solids  and  numerical 
integration.  The  insight  gained  in  this  education  has  been  valuable  and  will 
be  applicable  in  any  further  work  done  in  vibration  and  acoustics. 

In  addition  to  the  alternate  application  of  this  technique  presented  in 
the  Results  section,  there  is  potential  application  possibilities  in  the  related 
problems  of  debonding  and  seam  integrity. 

Debonding  occurs  when  the  coating  layer  does  not  properly  adhere 

to  the  steel  plate  substrate  layer.  This  can  be  a  disastrous  circumstance 

since  the  portion  of  the  polymer  layer  that  is  not  secure  to  the  hull  may 

separate  in  the  presence  of  viscous  flow  (as  the  boat  moves  through  the 

water)  and  collide  with  apparatus  aft  of  the  area  from  which  the  coating 

section  departed.  This  is  another  situation  which  would  not  become 

apparent  until  well  after  the  coating  application  has  been  completed.  The 

technique  presented  here  would  be  able  to  investigate  the  bond  integrity 

since  part  of  the  composite  signal  is  due  to  the  reflection  off  of  the 

interface  between  the  polymer  layer  and  the  steel  plate.  The  reflection 

coefficient,  which  described  the  percentage  of  the  energy  returned  into 

the  coating  after  reflection,  is  dependent  on  the  type  of  bond  at  the 

interface  and  is  opposite  in  sign  for  the  case  of  a  completely  debonded 

interface.  There  is  also  additional  avenues  of  pursuit  available  for  exploring 

this  problem  with  the  use  of  rotational  shear  (torsional)  waves.  Since  the 
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particle  motion  in  this  type  of  wave  is  parallel  to  the  interface  plane  there 
would  be  considerable  impact  on  the  reflection  coefficient  if  the  bond  was 
poor  or  non-existent. 

Seam  integrity  is  a  description  of  the  adherence  of  adjacent  coating 
strips  to  each  other,  (see  introduction  for  a  brief  description  of  the 
application  process).  Poor  seam  integrity  is  problematic  when  the  boat  is 
underway  since  it  may  lead  to  failures  that  cause  section  of  the  coating  to 
depart  from  the  hull  by  allowing  flow  forces  from  the  medium  to  peel  the 
edge  of  a  coating  section  up  and  eventually  force  an  entire  section  to  tear 
away  from  the  hull.  The  technique  presented  here  is  applicable  to 
investigating  this  failure  mode  since  part  or  all  of  the  composite  signal 
received  at  the  surface  of  the  coating  would  travel  across  a  seam  if  the 
transmit  and  receive  elements  were  positioned  on  opposite  sides  of  a 
seam.  Poor  seam  adherence  would  manifest  itself  primarily  by  reduction  in 
signal  amplitudes  due  to  the  reflection  of  energy  at  the  interface  between 
poorly  adhered  adjacent  coating  strips.  There  may  also  be  some  distortion 
in  the  case  of  a  marginal  bond  that  would  introduce  new  propagation 
dynamics  into  the  signals. 

These  alternate  applications  of  the  technique  may  be  a  valuable  use 
of  the  technology  developed  here  since  they  might  provide  a  binary 
("yes"  or  "no")  method  for  evaluating  a  given  instance  of  coating 
application. 

An  additional  consideration  for  the  application  of  this  research  to 
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usable  industry-grade  technology  is  the  work  environment.  The 
experiments  involved  in  this  research  effort  were  conducted  in  a 
laboratory  and  as  such  were  afforded  favorable  conditions  that  might  not 
exist  in  an  actual  environment,  e.g.  a  naval  construction  shipyard. 
Particularly,  dust,  noise,  and  vibration  are  variables  that  would  affect  the 
operability  of  such  a  system.  To  overcome  these  problems  would  involve 
significant  refinement  of  the  technology  involved  to  improve  its  robustness 
not  only  in  terms  of  physical  strength  but  operational  simplicity.  The 
current  state  of  the  technology  requires  a  knowledgeable  operator 
conducting  the  experiment  with  vigilance  which  would  not  likely  be 
available  in  the  shipyard  environment.  However,  it  might  be  feasible  to 
achieve  this  type  of  robustness  for  a  binary-decision  type  of  system 
described  above. 
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Appendix  A 


LDV  System  Details 

This  research  effort  used  a  Laser  Doppler  Velocimeter  already 
resident  at  Georgia  Tech.  It  used  a  35mW  Helium-Neon  laser  as  its  light 
source  and  optical  fibers  to  deliver  the  light  to  the  test  surface. 

The  basic  operational  principle  of  the  Laser  Doppler  Velocimeter  is  to 
detect  the  change  in  frequency  of  the  light  reflected  from  a  vibrating 
surface.  The  amplitude  and  phase  of  the  surface  velocity  is  determined 
from  the  measured  Doppler  frequency  shift.  These  measured  quantities 
can  then  be  processed  to  determine  the  particle  motion  of  the  surface 
which  is  being  illuminated  by  the  focused  laser  spot.^ 

The  general  layout  of  the  hardware  is  shown  in  figure  15.  In  the 
figure  the  components  are  labeled  by  letters  from  A-N.  A  description  of 
each  piece  and  its  function  is  given  here. 

A:  Helium-Neon  laser;  generates  a  coherent  beam  with  a  wavelength  of 

632.8  nanometers 

B:  Beamsplitter  cube;  reflects  approximately  4%  of  the  light  to  form  a 

reference  beam 

C:  Acousto-Optic  modulator;  a  quartz  crystal  driven  a  40MHz  oscillator 

that  is  adjusted  to  split  the  light  into  many  beams,  two  of  which  are 
approximately  equal  in  amplitude  but  differ  in  frequency  by  the  oscillator 

^  Hyun-Gwon  Kil,  Jacek  Jarzynski,  Yves  Berthelot.  "Wave  decomposition  of  the  vibrations  of  a 
cylindrical  shell  with  an  automated  scanning  laser  vibrometer".  Journal  of  the  Acoustical  Society  of 
America,  November  1998. 
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frequency. 

D:  mirror;  reflects  the  diverging  beams  coming  out  of  the  AOM  (C) 

E:  mirror;  reflects  only  one  beam 

F:  mirror;  reflects  only  one  beam,  shifted  by  40MHz  from  the  beam 

that  is  being  reflected  by  mirror  E 

G-I:  beam-fiber  coupler;  focuses  airborne  light  beam  into  the  core  of  an 
optical  fiber. 

G-  used  for  the  reference  leg  in  out  of  plane  detection 
H-  used  as  signal  leg  for  out  of  plane  detection  or  as  one  of  the  legs 
in  the  in-plane  detection 
I-  used  as  other  leg  in  the  in-plane  detection 
J :  beamsplitter  cube;  used  for  re-combination  of  the  reference  leg  and 

signal  leg  in  out  of  plane  detection 

K:  focusing  lens;  used  to  focus  the  recombined  interference  signal  onto 

the  active  surface  of  the  optical  detector  (avalanche 
photodiode) 

L:  Avalanche  photodiode;  used  to  convert  the  optical  signal  beam  into 

an  electrical  signal  for  digitization,  storage,  and  analysis 
M:  optical  probe  head;  used  to  focus  and  aim  the  light  beams  onto  the 

subject  test  surface  and  collect  the  reflected  (specular  and 
diffuse)  light  from  the  surface 

N :  test  surface;  the  object  whose  motion  is  to  be  detected 
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c::>  K 


L 


figure  18:  LDV  setup  schematic 


The  optical  probe  head  utilized  GRIN  lenses  mounted  on  adjustable 
springs.  The  optical  fibers  were  positioned  using  small  diameter  tubing  and 
custom-built  holders.  See  figure  19. 
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figure  19:  close-up  view  of  the  optical  probe  head 


For  more  details  on  the  operation  of  the  LDV  system  and  Laser 
Doppler  Vibrometry,  please  refer  to  "Wave  decomposition  of  the 
vibrations  of  a  cylindrical  shell  with  an  automated  scanning  laser 
vibrometer"  by  Hyim-Gwon  Kil,  Jacek  Jarzynski,  and  Yves  Berthelot. 
Also,  The  Laser  Doppler  Technique  by  L.E.  Drain  (New  York:  John  Wiley 
and  Sons,  1980)  is  an  excellent  source  for  information  on  Laser  Doppler 
Vibrometry, 
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Appendix  B 

Simplex  Computer  Code 


This  code  was  programmed  from  the  FORTRAN  code  listed  in  the 
Numerical  Recipes^  book,  based  on  the  work  of  Nelder  and  Mead.  A  test 
Simplex  solution  is  presented  here  for  a  2-dimensional  function.  A  plot  has 
been  made  of  the  solution  simplexes  to  aid  in  visualizing  the  search  path  of 
the  algorithm. 


‘  Numerical  Recipes:  The  Art  of  Scientific  Computing.  (New  York:  Cambridge  University  Press,  1986). 
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X  value 


figure  20:  sample  Simplex  result  for  a 
simple  function 
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function:  (x-1.4)  +(y-0.3) 


The  Simplex  solution  code  was  translated  from  the  FORTRAN  code 
in  Numerical  Recipes  and  is  presented  here. 

function  simplex2000(p,ftol4tmax,alpha,beta,gamma) 

%  From  Numerical  Recipes  page  292-  Nelder  and  Mead  alogrithm 
%  programmed  March  13,  2000 

% 

%  usage:  simplex2000(p,ftol,itmax,alpha,beta,gamma) 

%  set  p=0  for  automatic  generation 
%  omit  alpha,  beta,  gamma  for  default  values 

%  debugging  option-  commentary 
do_comments=0; 

%  code  for  visually  tracking  progression 

newFig=figure('unitsVnormVposition',[0,0,l,l]); 

set(gcf,'units','pixels'); 

figPos=get(gcf,'position'); 

scrnSize=figPos(3:4); 

xlabel('x  value');ylabel('y  value'); 

title('Simplex  Progression') 

grid 

%delete(newFig); 

%  Create  an  initial  simplex  if  one  was  not  passed  in 
if  p==0, 

p=[0.7,0.2;0.8,0.24;0.72,0.46]; 

end; 

%  Set  alpha,  beta,  and  gamma  if  not  passed  in 
if  nargin==3, 
alpha=l; 
beta=0.5; 
gamma=2.0; 

end; 

ndim=size(p,l)-l; 

mpts=ndim+l; 
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iter=0; 

for  index=l:size(p,l)/ 

y(index)=sfn([p(index,l),p(index,2)]); 

end; 

rtol=2’^ftol; 

while  abs(rtol)>ftol,  %  [10] 
ilo=l; 

if  y{l)>y(2), 
ihi=l; 
inhi=2; 

else 

ihi=2; 

inhi=l; 

end 

for  index=l:mpts,  %[11] 

if  y(index)<y(ilo),ilo=index;,end; 
if  y  (index)  >y(ihi), 
inhi=ihi; 
ihi=index; 

elseif  y  (index)  >y(inhi), 

if(index~=ihi),inhi=index;,end; 

end; 

end;  %[11] 

rtol=2*abs(y(ihi)-y(ilo))/  (abs(y(ihi))+abs(y(ilo))); 
if  iter==itmax 

dispCexceeding  max  iterations') 

P 

[sum(p(:,l)/3),sum(p(:,2)/ 3)] 
pause 

end; 

iter=iter+l; 

for  jndex=l:ndim,  %  [12] 

pbar(jndex)=0; 
end;  %  [12] 


50 


%  [14] 


for  index=l;mpts, 
if  index~=ihi, 

for  jndex=l:ndim,  %  [13] 

pbar(jndex)=pbar(jndex)+p(index,jndex); 
end;%  [13] 

end; 

end;%  [14] 

if  do_comments==l, 

dispCextrapolate  by  alpha  through  face')  %  — > 

end; 

for  jndex=l:ndim,  %  [15] 

pbar(jndex)=pbar(jndex)/  ndim; 
pr(jndex)=(l+alpha)*pbar(jndex)-alpha’^p(ihi,jndex); 
end;%  [15] 

if  do_comments==l, 

dispCevaluate  at  reflected  point')  %  --> 

end; 

ypr=sfn(pr); 
if  ypr<y(ilo), 

if  do_comments==l, 

dispCresult  better  than  best,  so  extrapolate  by  gamma') 

%  -> 

end; 

for  jndex=l:ndim,  %  [16] 

prr(jndex)=gamma’*^pr(jndex)+(l-gamma)'^pbar(jndex); 
end;  %  [16] 

yprr=sfn(prr); 
if  yprr<y(ilo), 

if  do_comments==l, 

dispCadditional  extrapolation  succeeded,  replacing 

highest  point')  %  --> 
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end; 


for  jndex=l:ndim,  %  [17] 

p(ihi,jndex)=prr(jndex) ; 
end;  %  [17] 

y(ihi)=yprr; 

else 

if  do_comments==l; 

dispCadditional  extrapolatio  failed,  still  using 

reflected  point')  %  --> 
end; 

for  jndex=l:ndim,  %  [18] 

p(ihi,jndex)=pr(jndex); 
end;  %  [18] 

y(ihi)=ypr; 

end; 

elseif  ypr>y(inhi), 

if  do_comments==l, 

dispCreflected  point  is  worse  thans  second  highest') 

%  -> 

end; 

ifypr<y(ihi), 

if  do_comments==l, 

dispCbetter  than  highest,  replacing')  %  — > 

end; 

for  jndex=l:ndim,  %  [19] 

p(ihi,jndex)=pr(jndex); 
end;  %  [19] 

y(ihi)=ypr; 

end 

for  jndex=l:ndim,  %  [21] 

prr(jndex)=beta’^p(ihi,jndex)+(l-beta)’^pbar(jndex); 
end;  %  [21] 
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yprr=sfn(prr); 
if  yprr<y(ihi). 


if  do_coinments==l, 

dispCcontraction  was  improvement')  %  — > 

end; 

for  jndex=l:ndim,  %  [22] 
p(ihi,jndex) =prr  (jndex) ; 
end;%  [22] 

y(ihi)=yprr; 


if  do_comments==l, 

dispCcontracting  around  lowest  point')  %  — > 

end; 


for  index=l:mpts,  %  [24] 

if  index<ilo, 

for  jndex=l:ndim,  %  [23] 

pr(jndex)=0.5*(p(index,jndex)+p(ilo,jndex)); 

p(index,jndex)=pr(jndex); 
end;  %  [23] 

y  (index)=sfn(pr) ; 

end; 

end;  %  [24] 

end; 

else 

for  jndex=l:ndim,  %  [25] 

p(ihi,jndex)=pr(jndex); 
end;  %  [25] 

y(ihi)=ypr; 

end; 


ifndim==2, 
hold  on 

fill(p(:,l),p(:,2),rand(l,3)); 

text((sum(p(:,l)))/3,(sum(p(:,2)))/3,num2str(iter)); 
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hold  off 


if  ndun==3, 
hold  on 

hold  off 


end;  %  [10] 

disp(['final  values  are: ']) 

P 

[sum(p(:,l )  /  3),sum(p(:,2)  /  3)] 

dispC  ');disp(['on  num2str(iter),'  iterations']) 

function  out=sfn(inputArgs) 
if  length(inputArgs)==2, 

x=input Args(l  );y =input  Args(2); 
out=(x-1.4)^2+(y-0.3)^2; 
elseif  length(inputArgs)==3, 

x=inputArgs(l);y=inputArgs(2);z=inputArgs(3); 

out=(x-1.4)^2+(y-0.3)^2+(z+2.7)^3; 

end; 
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Appendix  C 


Weaver /Hankel  Model  Code 


The  Hankel  Inverse  Transform  code  derived  from  Weaver's  paper 
presented  here. 


%  Calculates  propagation  of  elastic  waves  along  the  surface.MillerP.m 
%  with  z=0.  Calculates  the  generation  of  elastic  waves  by  a  normal 
%  load  on  the  surface  of  a  layer  on  rigid  back.  Based  on  the 
%  integral  transform  solution  of  Weaver. 

%  Includes  integration  using  Simpson  1/3  rule. 

%  Calculation  of  propagation  of  the  transient  signal 
%  Uses  accelerometer  signal  as  the  drive  -  pulse(n) 

clear  all; 

load  drive2p.dat; 

%  program  to  create  a  one  cycle  sine  wave  pulse, 

%  which  is  the  drive  signal. 

for  n= 1:256 
pulse(n)=0; 
end 

dtau=2.4e-5; 

for  n=  1:256 

pulse(n)=drive2p(n); 

end 

%create  time  array 

for  n=l:256 

time(n)=(n-l)’^dtau; 

end 

%Fourier  transform  of  the  drive  pulse 
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Fc=fft(pulse); 


%  create  propagator  array  for  given 
%  position  of  observation  point  and  drive  point 
a=0.0025; 
r-0.064; 

%  z=0.0; 

%  properties  of  layer  material 
L=0.05; 

%  frequency  loop 

delfr=l/(dtau*256); 

for  1=1:128 
uzarr(l)=0; 
end 

for  1=2:128 
sho=l 

frarr(l)=(l-l)’*^delfr; 

freq=frarr(l); 

clr=2.418e+l*log(freq)  +  7.303e+l; 
cli=8.62inog(freq)  -  2.152e+l; 
csr=9.535*log(freq)  +  2.857e+l; 
csi=3.332’^log(freq)  -  8.019; 
cl=dr+i*cli; 
cs=csr+i’^csi; 

omeg=2’^pi’*'freq; 

kl=omeg/cl; 

k2=omeg/cs; 

rk2=real(k2); 

k2p=k2/rk2; 

klp=kl/rk2; 

delk2p=l/500; 

ap=a’^rk2; 

rp=r’^rk2; 
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Lp=L=^rk2; 

%  integrand  array  for  uz 
clear  zetapt  vlpt  besst  uzintt 


zmax=round(61-0.3*l); 

zetapt=0:delk2p:zmax; 

max=zmax/delk2p+l; 

uzintt=zeros(l,max); 

vlpt=sqrt(klp*klp-zetapt.*zetapt); 

besst=bessel(l,ap’^zetapt).’^bessel(0,rp’^zetapt); 

%check  maximum  value  of  vlpt 

if  -imag(vlpt(max))*Lp<150 
ind=l; 
else 
ind=2; 
end 

if  md==2 

mdiv=round(150/ (Lp*delk2p)); 
end 

%  develop  full  or  asymptotic  expression  for  the 
%  z  and  r  integrands.If  ind=l  use  full  expression 
%  for  uz  and  ur  over  the  entire  range  of  zetapt. 

%  If  ind=2  use  full  expression  for  zetap  less  than 
%  or  equal  to  mdiv  and  the  asymptotic  expression 
%  for  zetap  greater  than  mdiv. 

if  md==l 

fzetat=(2’^(zetapt.*zetapt)-k2p*k2p); 
vlt=sqrt(zetapt.*zetapt-klp*klp); 
v2t=sqrt(zetapt.*zetapt-k2p*k2p); 
v2pt=sqrt(k2p’^k2p-zetapt.*zetapt); 
clpt=2*vlt.  /fzetat; 

c3pt=clpt.’^((zetapt.*zetapt).*cos(Lp*v2pt))-vlt.*cos(LpVlpt); 
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dlpt=(clpt.’^(zetapt.’^zetapt).*(cos(Lp'^vlpt).*cos(Lp’^v2pt))); 
dlpt=dlpt+(clpt.*(vlt  *v2t)).*(sin(LpVlpt).*sin(Lp*v2pt))-vlt; 
d2pt=(clpt.’^v2t).*(zetapt.’^zetapt)-(vlt  *v2t).*(cos(Lp’^v2pt).’^cos(Lp'^vlpt)); 
d2pt=d2pt-(zetapt.*zetapt).’^(sin(LpV2pt).’^sin(LpVlpt)); 

uznumt=((zetapt  *zetapt).*(dlpt.'^sin(LpV2pt))-vlt.’^(d2pt.’^sin(Lp^vlpt))) 

.’^(i’*^((zetapt.*zetapt).*clpt-vlt)); 

uzdent=c3pt.*(fzetat.’^d2pt-2*(zetapt.’*'zetapt).’^(v2t.*dlpt)); 

%  integrand  arrays  for  uz  and  ur 

uzintt=(uznumt./  uzdent).*besst; 
clear  fzetat 

clear  vlt  v2t  v2pt  dpt  c3pt  dipt 
clear  d2pt  uznumt  uzdent 
end 

if  md==2 

zetap=zetapt(l:mdiv); 

zetapa=zetapt(mdiv+l:max); 

bess=besst(l  :mdiv); 
bessa=besst(mdiv+l:max); 

fzeta=(2*(zetap.*zetap)-k2p*k2p); 

vl =sqrt(zetap .  ’^zetap-klp’^klp); 

v2=sqrt(zetap.*zetap-k2p*k2p); 

vlp=sqrt(klp*klp-zetap.*zetap); 

v2p=sqrt(k2p’^k2p-zetap.’^zetap); 

clp=2’^vl./fzeta; 

c3p=clp.*((zetap.*zetap).’^cos(LpV2p))-vl.*cos(LpVlp); 

dlp=(clp.“^(zetap.*zetap).’^(cos(Lp*vlp).*cos(Lp’^v2p))); 

dlp=dlp+(clp.*(vl.*v2)).*(sin(Lp’^vlp).*sin(Lp’^v2p))-vl; 

d2p=(clp.’*'v2).*(zetap.*zetap)-(vl.*v2).*(cos(LpV2p).*cos(LpVlp)); 

d2p=ci2p-(zetap.*zetap).*(sin(LpV2p).’^sin(Lp*vlp)); 

uznum=i’^((zetap.*zetap).*(dlp  *sin(Lp*v2p))-vl  *(d2p.“^sin(Lp'^vlp))) ... 
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.“^((zetap.*zetap).’^clp-vl); 

uzden=c3p.’'(fzeta.’'d2p-2’'(zetap.’'zetap).’'(v2.^dlp)); 

%  integrand  arrays  for  uz  in  the  full  region 

uzinttl=(uznum./uzden).*bess; 
uzintt(l  :mdiv)=uzinttl ; 

clear  zetap  fzeta  bess 

clear  vl  v2  vlp  v2p  clp  c3p  dip 

clear  d2p  uznum  uzden  uzinttl 

fzetaa=(2*(zetapa.*zetapa)-k2p*k2p); 

vla=sqrt(zetapa.*^zetapa-klp’^klp); 

v2a=sqrt(zetapa.’^zetapa-k2p’*^k2p); 

vlpa=sqrt(klp*klp-zetapa.’^zetapa); 

v2pa=sqrt(k2p*k2p-zetapa.*zetapa); 

clpa=2’^vla./fzetaa; 

uznuma=((zetapa.’'zetapa).’'clpa-vla).’'((zetapa*zetapa).’'clpa  ... 

.=^exp((abs(vlpa)-abs(v2pa))*Lp)  -  via); 
uzdena=((zetapa.’'zetapa).’^clpa.’'exp((abs(vlpa)-abs(v2pa))='Lp)-vla)  ... 

.’^(fzetaa-2V2a.*(zetapa.*zetapa).*clpa); 

uzintta=(uznuma./uzdena).’^bessa; 

uzintt(mdiv+l  :max)=uzintta; 

clear  zetapa  fzetaa  bessa  bessra 
clear  via  v2a  vlpa  v2pa  clpa  c3pa  dlpa 
clear  d2pa  uznuma  uzdena  uzintta 

end 

%  integrate  the  arrays  uzintt 
n=zmax’^500; 

uzfr=uzintt(l)+4’'sum(uzintt(2:2:n))+uzintt(n+l); 

uzfr=uzfr+2’^sum(uzintt(3:2:n)); 

uzfr=(delk2p/3)’^uzfr; 
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uzarr(l)=uzfr; 

vzarr(l)=i*omeg*uzarr(l); 

end 

frarr(l)=0; 
uzarr(l  )=uzarr(2) ; 
vzarr(l)=vzarr(2); 

%  create  inverse  Fourier  transform  of  ur,uz  and  vr,vz  signals 
for  n=l:128 

uzinvs(n)=uzarr(n)’^Fc(n); 

vzinvs(n)=vzarr(n)*Fc(n); 

end 

for  n=130:256 
k=n-130; 

uzinvs(n)=conj(uzmvs(128-k)); 

vzinvs(n)=conj(vzinvs(128-k)); 

end 

uzinvs(129)=uzinvs(l); 

vzinvs(129)=vzinvs(l); 

%Take  inverse  Fourier  transform  to  generate  uz  and  vz, 
%and  ur  and  vr  signals  at  r,z=0 

siguz=if ft(uzmvs) ; 
sigvz=ifft(vzinvs); 

%plot(time,.01’^pulse/b') 

%hold 

%plot(time,real(siguz),'r') 

plot(time(l:100),real(-sigvz(l:100))/r') 

end 


60 


Appendix  D 


Finite  Element  Modelling 

In  section  VII,  Comments  and  Conclusions,  some  discussion  and 
results  were  presented  for  work  done  on  cavity  detection  capability  with 
the  current  system.  In  addition  to  this  experimental  work  there  was  some 
finite  element  modelling  to  investigate  the  effects  and  detectability  of  an 
occlusion  in  a  sample  of  the  type  analyzed  in  the  present  work.^  The 
modelling  was  intended  to  analyze  the  influence  of  a  cavity  on  the  motion 
at  the  interface  between  the  polymer  material  and  medium,  in  the  case  air 
(free  surface  of  sample). 


figure  21:  local  impedance  discontinuity  causes  flow  towards  cavity 


^  Caspall,  Jayme.  "Dynamics  of  Cavities  in  Layered  Media."  Presentation  given  as  part  of  the  In-Situ 
Coating  Program  review  presented  on  19  August  1997. 
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Figure  21  depicts  the  way  that  mass  in  the  sample  will  tend  to  flow 
towards  the  cavity  in  the  layer.  Since  there  is  a  mass  flow  induced  by  the 
impedance  discontinuity  there  will  be  an  associated  resonance  from  the 
effective  spring  caused  by  the  impedance  and  mass  flow.  A  resonance  of 
this  type  is  of  interest  in  the  case  where  sensors  are  attached  to  the  free 
surface  of  the  sample  layer.  When  positioned  as  a  linear  array,  the  sensors' 
output  will  be  affected  by  the  noise  generated  from  the  nonuniform 
surface  motion. 


incident  plane  wave 


sensor  output 


00000°  °oooooooooooo 


sensor  position 


figure  22:  nonuniform  sensor  output 


Finite  element  models  were  created  using  Sara  2D® .  One  of  them 
was  used  to  investigate  the  effect  of  cavity  depth  on  resonance  and  the 
other  the  effect  of  cavity  height  on  resonance.  In  the  former,  a  cavity  is 

"  "Sara  2D"  by  BBN  Systems  and  Technologies,  Union  Station,  New  London,  CT  06320-6147. 

62 


inserted  at  the  free  surface  and  its  size  in  the  axial  direction  varied.  In  the 
latter  the  cavity  is  inserted  at  the  constrained  surface  (polymer  lab  bonded 
to  steel  plate)  and  grows  axially. 


REZCON  Resonance  vs.  CavityDepth 


0  0.05  0.1  0.15  0.2  0.25  0.3  0.35  0.4 


depth  of  cavjty(in.) 

figure  23:  effect  of  cavity  depth  on  resonance;  vertical  axis 
shows  reduction  in  resonance  frequency  compared  to  no  cavity,  % 
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REZCON  Resonance  vs.  CavityHeighl 


height  of  bubble  (in.) 

figure  24:  effect  of  cavity  height  on  resonance;  vertical  axis 
shows  reduction  in  resonance  frequency  compared  to  no  cavity,  % 

A  third  FEM  was  used  to  see  the  effect  of  cavity  radius  on 
resonance. 
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R  EZCON  R  es  onanc  e  vs .  Cavity  Radi  us 


figure  25:  effect  of  cavity  radius  on  resonance;  vertical  axis 
shows  reduction  in  resonance  frequency  compared  to  no  cavity,  % 

The  finite  element  modelling  was  used  for  detectability  estimation 
using  a  two-dimensional  plane  strain  model.  The  load  was  a  unit  line  force 
and  the  cavity  was  treated  as  an  infinite  cylindrical  free  boxindary  within 
the  coating.  The  velocity  was  calculated  for  various  distances  from  the  line 
force  and  covered  a  range  commensurate  with  actual  measurements.  Plots 
are  shown  comparing  the  "with-cavity"  case  to  the  "without-cavity"  case 
and  the  calculations  are  repeated  for  various  cavity  diameters  and  various 
cavity  depths. 
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figure  26:  FEM  used  for  detectability  estimation 


Norm.  surf.  vel.  forf  =  2000  Hz 


distanoe  Irom  shaker  (in.) 


Difference  Between  With  and  Without  Cavity 


distance  from  shaker  (jn.) 


figure  27:  use  of  FEM  results 
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The  model  results  are  shown  in  figure  17.  At  a  given  frequency,  the 
peak  relative  difference  in  the  amplitude  of  the  normal  surface  velocity 
between  the  with-cavity  and  the  without-cavity  case  is  calculated  for  each 
cavity  diameter  and  depth.  The  calculated  difference  is  compared  to  he 
experimentally  determined  minimum  detectable  difference. 

Figure  28,29,  and  30  are  cavity  detectability  matrices  from  model 
calculations.  Each  shows  the  detectability  status  for  combinations  of 
frequency,  cavity  diameter,  and  cavity  depth.  Figure  28  shows  results 
based  on  current  capabilities,  figure  29  shows  results  with  an  attainably 
better  system,  and  figure  30  shows  a  cumulative  report. 
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figure  28:  cavity  detectability 
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figure  29:  cavity  detection  based  on  improved  system 
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1/16"®  ©©©©©© 

0.03"  055“0.50"0.75"  1.00“  155"  1.5" 
depth 


High  Resolution  Goal 
Minimum  Detectable  Signal  2% 

1A"®  ®  ®  ®  ®  ®  ® 

|/t6"®  •  ®  •  ®  ®  ® 

I  1^®  ®  ®  ®  ®  ®  ® 

1/16“®  ®  @  ©  ©  ©  © 

0.03"0.25"0.50"0.75"  1.00"  1.25"  1.5" 
depth 


’^Note:  4  kHz  is  the  highest  frequency  allowed  by  the  model.  The 
fundamental  cavity  resonance  is  higher  than  4  kHz.  The  1/4”  cavity 
has  the  iowest  resonance  frequency  of  all  the  cavities  in  the  study  at 
approximately  7  kHz.  Scattering  results  for  the  case  of  a  spherical 
cavity  in  an  infinite  homogeneous  medium  suggests  that  the  resonance 
provides  the  largest  signals  for  detection. 

figure  30:  cumulative  detectability 
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Some  analytical  calculations  were  done  to  look  at  the  effects  of 
detecting  scattered  signals  from  an  occlusion  at  resonance  versus  not  at 
resonance.  The  scattering  calculations  are  based  on  the  formulation  of 
Graff  with  plane  wave  incidence.^  In  reference  to  the  plots,  the  cavity  is 
centered  at  x=0  and  the  detection  threshold  is  provided  as  a  reference 
only.  The  scattering  formulation  assumes  an  infinite  medium  and  cannot 
provide  accurate  detection  estimates  for  the  case  of  a  coating  layer  on  a 
steel  plate. 


Scattering  Calculations  On  Resoncmce  and  Off 


Scattered  Displacement  On  a  Plane 

Located  1"  from  0.25"  Diameter  Cavity 
Normal  Component 


0  os  1  IS  2  2S  3  3S  4  4S  S 

X  Onches  flom  cavity  axis) 


figure  31:  scattering  calculations  on  resonance  and  off 


’  Graff,  Karl  F.  Wave  Motion  in  Elastic  Solids.  (New  York:  Dover  Publications,1973). 
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